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ABSTRACT 

We find a remarkably enhanced production rate in star clusters (relative to the 
field) of very short period, massive double-white-dwarf stars and of giant-white 
dwarf binaries. These results are based on A-body simulations performed with 
the new GRAPE-6 special purpose hardware and are important in identifying and 
characterizing the progenitors of type la supernovae. The high incidence of very 
close double-white-dwarf systems is the result of dynamical encounters between 
(mostly) primordial binaries and other cluster stars. Orbital hardening rapidly 
drives these degenerate binaries to periods under ~ 10 hours. Gravitational 
radiation emission and mergers producing supra- Chandrasekhar objects follow 
in less than a Hubble time. If most stars are born in clusters then estimates 
of the double white dwarf merger rates in galaxies (due to cluster dynamical 
interaction) must be increased more than tenfold. A majority of the Roche lobe 
overflow giant-white dwarf binaries are not primordial; they are produced in 
exchange reactions. Most cases resulted in a common-envelope and formation of 
a double-white-dwarf binary rather than Supersoft X-ray sources leading possibly 
to a type la supernova. 

Subject headings: stellar dynamics — methods: N-body simulations — type la 
supernovae — globular clusters: general — open clusters and associations: general 



1. Introduction 



Type la supernovae (SNela) have recently been used to demonstrate that the Universe 
is apparently not only expanding, but also accelerating (Riess et al. 1998; Perlmutter et 
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al. 1999). If this remarkable discovery survives careful scrutiny it has profound implications 
for cosmology physics (see Leibundgut 2001 for a review). However, before accepting such 
a revolutionary change in our view of the universe, it is critical to investigate all challenges 
to the acceleration interpretation of the data. A thorough compilation of those challenges is 
given by Riess (2000). 

Riess (2000) cites evolution, dust, gravitational lensing, measurement biases, selection 
biases and alternative cosmological models as potential challenges to the acceleration inter- 
pretation of the Type la supernova (SNIa) observations. He concludes that the primary 
source of reasonable doubt is evolution (could SNela at redshift z = 0.5 be intrinsically 
fainter than nearby SNela by 25% ?). We simply don't know for certain what kind of star 
(or stars?) give rise to these explosions, and hence whether to expect systematic variations 
in SNIa luminosities with z. 

Not knowing the progenitors of SNela is embarrassing because of the significant empirical 
corrections one must apply to supernova luminosities, based on their light curves, to get 
distances (Phillips 1993). If one could determine the SNIa progenitor with a high degree 
of reliability, one could build increasingly sophisticated supernova light curve models. This 
would lead to a fundamental understanding of the physical behavior of one of the most 
crucial standard candles in cosmology. 

Two competing SNIa models exist: merging double-white-dwarfs (DWDs) and accreting 
single white dwarfs (ASDs) in close binary systems (see Yungelson & Livio 2000 for a review, 
including the pros and cons of each flavor of each model). In both cases the model involves the 
thermonuclear disruption of a white dwarf, most likely of carbon-oxygen (CO) composition, 
when its mass reaches, or exceeds, the critical Chandrasekhar mass. However, Saio & Nomoto 
(1998) have raised the possibility that supra- Chandrasekhar mass-accreting white dwarfs will 
undergo accretion-induced collapse (AIC) and the formation of a neutron star, rather than 
a SNIa. An important example of why it is so critical to know the progenitors of SNela 
is the following. ASDs accreting helium result in the accumulation of a He layer and an 
edge lit detonation; the peak luminosity-light curve shape relation of such objects may vary 
significantly with metallicity and hence redshift (Tout et al. 2001). 

It is clearly important to calculate the expected incidence of DWDs and ASDs in the 
stellar populations of the different types of galaxies where SNIa occur. This has been done 
for field stars (Yungelson et al. 1996, for example), but hardly for the environs of clusters. 
We have begun this study, and find that dynamics dramatically alter binary populations and 
characteristics, including those of DWDs and ASDs. Such an effect has been predicted in 
the past (Chen & Leonard 1993) but has yet to be tested by direct means. 
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Many, and possibly all, stars were born in a star cluster (Kraft 1983; Lada, Strom 
& Myers 1993), or at the very least a loose association (Makarov & Fabricius 2001, for 
example). Living in an environment with stellar densities of 10 2 stars pc -3 to as high as 
10 7 stars pc -3 can dramatically affect the evolution of stars. Physical collisions, and, in the 
case of binaries, exchange interactions or disruptions of orbits can radically alter the fates 
of cluster stars. The most direct way to model the evolution of a star cluster is with an N- 
body code in which the individual orbits of each star are followed in detail and the internal 
evolution of each star is also taken into account (Hurley et al. 2001). 

We have begun a study of the behavior of populous star clusters using a state-of-the-art 
A-body code in conjunction with the powerful GRAPE-6 special purpose computer (Makino 
2001). This will ultimately involve a large number of A-body simulations covering a wide 
range of initial conditions, e.g. metallicity, binary fraction, and stellar number density. A 
related study investigating the fate of planetary systems in star clusters is also ongoing 
(Hurley & Shara 2002). The remarkable evolution of close binaries comprised of one or two 
degenerate members is the focus of this early work. The increased importance of SNIa for 
cosmology makes these initial results particularly relevant to the scientific community. 

Our simulation method is presented in Section 2, and the double white dwarfs are 
described in Section 3. The accreting white dwarfs are detailed in Section 4. We briefly 
summarize our results in Section 5. 



2. Simulation Method 

We have carried out simulations with 22 000 stars and a 10% primordial binary frac- 
tion using the Aarseth NB0DY4 code (Aarseth 1999; Hurley et al. 2001). A GRAPE-6 
board located at the American Museum of Natural History hosts the code. This special pur- 
pose computer acts as a Newtonian force accelerator for A-body calculations, performing at 
0.5 Tflops for a 16-chip board (~ 30 Gflops per chip). 

The prescription used for single star evolution in NB0DY4 is described in Hurley, Pols & 
Tout (2000). A feature of this algorithm is the inclusion of metallicity as a free parameter, 
making it applicable to modelling clusters of all ages. It covers all stages of the evolution 
from the zero-age main-sequence (ZAMS) up to, and including, remnant stages such as the 
white dwarf cooling track. In terms of examining DWDs as SNIa candidates it is important 
to note that the WD initial-final mass relation found by Hurley, Pols & Tout (2000, see 
their Fig. 18) is well matched to observations. 

All aspects of standard binary evolution, i.e. non-perturbed orbits, are treated ac- 
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cording to the prescription described in Hurley, Tout & Pols (2002). This includes tidal 
circularization and synchronization of the orbit, mass transfer, and angular momentum loss 
mechanisms such as magnetic braking and gravitational radiation. The default input pa- 
rameters to this algorithm (listed in Table 3 of Hurley, Tout & Pols 2002) are adopted. In 
particular, this means that the common-envelope (CE) efficiency parameter is taken to be 
«ce = 3.0 which makes the outcome of common-envelope evolution similar to the alternative, 
and commonly used, scenario described by Iben & Livio (1993) with q;ce = 1-0. Common- 
envelope evolution occurs when mass transfer develops on a dynamical timescale. In the 
evolution algorithm this equates to the donor star, generally a giant, having an appreciable 
convective envelope and a mass-ratio, q, exceeding some critical value, g crit — 0.7. If the 
conditions for dynamical mass-transfer are met then the envelope of the giant overfills the 
Roche-lobes of both stars leaving the giant core and the secondary star contained within a 
common-envelope. Owing to orbital friction these will spiral together and transfer energy 
to the envelope with an efficiency ckce- If this process releases sufficient energy to drive off 
the entire envelope the outcome will be a close binary consisting of the giant core and the 
secondary, otherwise it leads to coalescence of the two objects. 

Various models for binary evolution exist (Portegies Zwart & Verbunt 1996; Tutukov & 
Yungelson 1996; Tout et al. 1997, for example) with each having the similar goal of provid- 
ing a sufficiently detailed description of binary behaviour while remaining computationally 
efficient. The main difference between these and the Hurley, Tout & Pols (2002) model 
referred to in this work, and incorporated in to NB0DY4, is the much improved treatment 
of tidal interactions present in the latter. Tidal friction arising from convective, radiative 
or degenerate damping mechanisms is modelled and necessarily the stellar spins, which are 
subject to tidal circularization and synchronization, are followed for each star. Also, by 
using the single star evolution algorithm of Hurley, Pols & Tout (2000) the Hurley, Tout & 
Pols (2002) binary prescription not only allows for a wider range of evolution phases, with 
many of these modelled in more detail and based on updated stellar models, but can also 
be used to evolve binaries of any metallicity. An additional difference is that models such 
as those of Portegies Zwart & Verbunt (1996) and Tutukov & Yungelson (1996) follow the 
Iben & Livio (1993) common-envelope scenario while Hurley, Tout & Pols (2002) utilise 
the scenario first described in Tout et al. (1997). Subtle variations also exist from model 
to model in the way that various aspects of mass transfer are dealt with. 

Gravitational radiation is an important process in the evolution of close binary systems 
because it provides a mechanism for removing angular momentum and driving the system 
towards a mass transfer state, possibly followed by coalescence. In the Hurley, Tout & 
Pols (2002) binary model orbital changes due to gravitational radiation are calculated using 
expressions based on the weak-field approximation of general relativity (Eggleton 2002) and 
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assuming the stars are point masses. This includes a strong dependence on the eccentricity 
of the orbit, although DWD systems are generally circular. 

In the dense environment of a star cluster it is possible for the orbital parameters of a 
binary to be significantly perturbed owing to close encounters with nearby stars. It is even 
possible for the orbit to become chaotic as a result of such an interaction. Modelling of such 
events has been considered in detail by Mardling & Aarseth (2001) whose work is included 
in NBODY4. Three-body and higher-order subsystems are also followed in detail (Aarseth 
1999, and references within). 

We include the outcome of four simulations in the results described below and in Tables 
1-3. The first two simulations were carried out assuming a metallicity of Z = 0.004, while 
the third and fourth simulations had Z = 0.02. In all other respects the initial conditions 
were identical for each simulation. Masses for single stars are chosen from the initial mass 
function (IMF) of Kroupa, Tout & Gilmore (1993) within the limits of 0.1 — 50M Q . For 
primordial binaries the total mass of the binary is chosen from the IMF of Kroupa, Tout & 
Gilmore (1991), as this was not corrected for the effect of binaries, between the limits of 
0.2 — 1OOM . The component masses are then assigned according to a uniform distribution of 
mass-ratio, taking care to ensure that the single star mass limits are not violated. Following 
Eggleton, Fitchett & Tout (1989) we take the distribution of orbital separations for the 
primordial binaries to be 



This distribution is symmetric in log a about a peak at a m and ranges from a minimum 
separation of (a m to a maximum of a m /(. We choose the constants ( and f3 to be 10~ 3 and 
0.33 respectively, and take a m ~ 30 AU, i.e. each separation, a, is within the range ~ 6 R Q to 
30 000 AU. The eccentricity of each binary orbit is taken from a thermal distribution (Heggie 
1975). Initial positions and velocities of the stars are assigned according to a Plummer model 
(Aarseth, Henon & Wielen 1974) in virial equilibrium. 

Each simulation started with 18 000 single stars and 2 000 binaries and was evolved to 
an age of 4.5 Gyr when ~ 25% of the initial cluster mass remained and the binary fraction 
was still close to 10%. In real time each individual simulation took approximately five 
days to complete, which corresponds to ~ 10 3 cluster crossing-times, and ~ 10 17 floating- 
point operations on the GRAPE board. Clusters are evolved subject to a standard three- 
dimensional Galactic tidal field (see Hurley et al. 2001) of standard, i.e. local, strength 




(1) 



where W G [—1,1], and uniformly distributed and k satisfies 



(2) 
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(Chernoff & Weinberg 1990, for example) so, in addition to mass loss from stellar evolution, 
mass is removed from the cluster when stars cross the tidal boundary and are lost to the 
Galaxy. Because the orbits of some bound stars may momentarily take them outside of the 
cluster, stars are not actually removed from the simulation until their distance from the 
cluster centre exceeds twice the tidal radius. At any moment in time less than 2% of the 
mass in a simulation is found to lie outside of the tidal radius and this has little effect on 
the evolution of the model (Giersz & Heggie 1997). Typically the velocity dispersion of the 
stars in these model clusters was 2kms~ 1 with a core density of 10 3 stars pc~ 3 . The density 
of stars at the half-mass radius is generally a factor of 10 less than this. These simulations 
are clearly in the open cluster regime. 



3. Double White Dwarfs 



Double white dwarf systems must form in all stellar populations with binaries. To 
form short-period systems isolated binaries must undergo considerable common envelope 
evolution to shed angular momentum and bring their degenerate components close together. 
Only when the orbital period of a DWD is less than about 10 hours will gravitational 
radiation force the system to merge in less than the age of the universe. A key result of 
this paper is that DWDs formed in clusters can be significantly hardened by interactions 
with passing stars. Alternatively, hardening of binaries before they reach the DWD stage, 
and/ or exchange interactions, can produce a short-period DWD that would not have formed 
in isolation. This greatly increases the number of DWDs which merge in less than a Hubble 
time. Heggie (1975) defined hard binaries to be sufficiently close that their binding energy 
exceeds the mean kinetic energy of the cluster stars. A binary is said to harden when an 
interaction with a third body removes energy from the orbit and thus reduces the separation. 
Three-body interactions may also lead to an exchange in which one of the binary components 
is displaced by the incoming third star. If an exchange does occur then the expelled star, 
generally the least massive, invariably leaves the three-body system altogether. Note that 
four-body interactions involving a second binary are also possible. The likelihood of a binary 
being the target of an exchange interaction scales linearly with the orbital separation of 
the binary (Heggie, Hut & McMillan 1996) and is also more likely to occur in the core 
of a cluster. Table 1 lists the characteristics of the merging DWDs formed during the four 
cluster simulations: epoch of formation, types of white dwarf, component and system masses, 
orbital period at formation epoch, gravitational inspiral timescale, and whether the binary 
is primordial or due to an exchange. All of these objects will merge in less than a Hubble 
time, creating a supra-Chandrasekhar object which may yield a SNIa. 
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From binary population synthesis alone we expect about 10 DWDs to be present 
amongst 2000 binaries of solar metallicity after 4 Gyr of evolution. This expected num- 
ber rises to 15 DWDs for Z = 0.004 which is mainly a reflection of the accelerated evolution 
timescale for low metallicity stars with masses less than ~ 9M : the main-sequence (MS) 
lifetime for stars of Population II composition is roughly a factor of two shorter than for 
stars of solar composition (see Fig. 4 of Hurley, Pols & Tout 2000). Of these DWDs, only 
~ 0.2 — 0.3 are expected to be a "loaded gun", i.e. to have a combined mass, M b , greater 
than the Chandrasekhar mass, M Ch ~ 1.44, and a merger timescale less than the age of the 
Galaxy (Hurley, Tout & Pols 2002). These population synthesis results refer to isolated 
binaries with initial conditions drawn from the same distributions as used in the iV-body 
simulations presented here. From Table 1 we see that there are, on average, four "loaded 
guns" per 2000 binaries, which is ~ 15 times more than expected in a binary population 
unaffected by dynamical interactions. We note CO-CO DWDs dominate, but that CO-ONe 
binaries are also common (5 of 16 binaries). The system total masses range from 1.49M Q to 
1.96M . This range in masses and compositions suggests diversity in the light curves and 
spectra of the merging objects, whether or not they make SNIa. 

Figure 1 highlights the evolution of two of the SNIa candidates formed in the simulations. 
The binary shown in the top panel of the figure is the DWD formed at 334 Myr in the first 
Z = 0.004 simulation (fourth entry in Table 1). This began as a primordial binary with 
component masses of 5.82M and 3.13M , an eccentricity of 0.74, and an orbital period of 
2138d. It underwent two common-envelope events, the first at 76 Myr when the 5.82M 
star filled its Roche-lobe on the asymptotic giant branch (by which time tides raised on the 
convective envelope of the primary had circularized the orbit) and became an ONe WD, and 
the second at 232 Myr when the 3.13M star filled its Roche-lobe on the first giant branch 
and became a naked helium star. The helium star subsequently evolved to become a CO WD 
at 299 Myr, losing some mass in the process which explains the decrease in binding energy. 
The DWD resided in the core of the cluster and experienced a strong perturbation to its 
orbit shortly after formation which caused the orbit to harden. The perturbation did not 
induce any noticeable eccentricity in to this already circular orbit. Gravitational radiation 
then removed orbital angular momentum from the system, causing the two WDs to spiral 
together, until they merged at ~ 630 Myr. Without the perturbation to its orbit the period 
of this DWD binary would not have been short enough for gravitational radiation to be 
efficient. 

The binary shown in the lower panel of Figure 1 is the DWD formed at 1 192 Myr in 
the second of the Z = 0.02 simulations (last entry in Table 1). This binary is the result 
of a 3-body exchange interaction which occurred at ~ 620 Myr. The binary involved in the 
exchange is a primordial binary that had experienced two common-envelope events and was 
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itself a DWD at the time of exchange (see the second entry in Table 2 for Simulation 4). 
The incoming star had a mass of 2.O2M and after spending a short time as a quasi-stable 
3-body system the least massive star, the O.82M WD in the primordial, was ejected from 
the system to leave a MS-WD binary with an eccentricity of 0.63 and an orbital period of 
14 125 d. This binary was significantly perturbed at 914 Myr (increasing the eccentricity to 
0.94) and then experienced a series of common-envelope phases (with the orbital eccentricity 
having been removed by tidal friction prior to the first CE event) before forming the SNIa 
candidate of interest at 1 192 Myr. It bears repeating that possible SNIa systems of this sort 
never occur naturally in the field. 

In Table 2 we list all DWDs which are not SNIa candidates. This is because the total 
system masses are less than M Ch , or the merger timescale is longer than a Hubble time, or 
both. The large numbers of objects demonstrate the rich variety of DWDs created in clusters. 
It is particularly remarkable that 93 of the 135 binaries listed in Table 2 were created in 
exchange interactions. This highlights a second key result of our paper: most of the DWDs 
in star clusters, and possibly those in the field if they were born in clusters, have progenitors 
with companions different from the ones we observe today. Furthermore, the orbital period 
and mass ratio distributions of such objects cannot be reliably predicted without including the 
effects of dynamics and thus, exchanges. This is particularly important in comparing the 
results of observational searches for "loaded guns" with theoretical population predictions 
(Saffer, Livio & Yungelson 1998). The number of DWDs formed from primordial binaries 
in the simulations matches well the predicted numbers from the non-dynamical population 
synthesis. This is also true of the SNIa candidates in Table 1 where the population synthesis 
predicted either or 1 candidate per simulation. 

The histogram of all DWD masses is shown in Figure 2. The distribution of DWDs which 
will merge in less than a Hubble time is quite similar to that of the DWDs with M b > M Ch . 
This is also clear from the histogram of all DWD periods, shown in Figure 3. These confirm 
that the hardening process does not preferentially act on more massive binaries. 

The color-magnitude diagrams of the clusters we have simulated are shown in Figures 4- 
5 for the Z = 0.004 and Z = 0.02 models at 1 Gyr intervals. In addition to the single star and 
binary main- sequences, several blue stragglers and cataclysmic variables (CVs, below and to 
the left of the MS) are visible in each simulation. The latter have hardly been searched for 
in open clusters, and this work suggests that systematic searches might be rewarded. Most 
remarkable are the DWDs, seen in profusion above the white dwarf cooling sequence, and 
the "loaded guns" shown as circled diamonds in each of the figures. 

Figure 6 illustrates the spatial distribution within the cluster at several epochs for the 
single stars, binaries and DWDs. This clearly demonstrates that the quest for equipartition 
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of energy is dominating the dynamical evolution: heavy objects segregate towards the inner 
regions of the cluster while lighter objects move outwards. Binaries are on average more 
massive than single stars and are therefore more centrally concentrated. Furthermore, the 
progenitor of any DWD was originally at least twice as massive as the cluster MS turn-off 
mass at the time the DWD formed, i.e. the DWDs form from the high-end of the binary 
IMF. This explains why DWDs are more centrally concentrated than the other binaries and 
means that they are most likely to be found in the core of the cluster, especially those with 
M b > M Ch . 



4. Accreting Single White Dwarfs 

Table 3 lists the possible ASDs generated in our simulations. As with the DWDs, we 
see that about 2/3 of these binary systems have been involved in exchange interactions. 
Most are giant-CO-WD pairs, but examples of He and ONe white dwarfs orbiting giants are 
also present. The number of ASDs is comparable to that predicted by population synthesis 
without dynamical interactions (~ 6 per simulation) but considering that only 1/3 of the 
ASDs evolved from primordial binaries it appears that dynamical interactions are destroying 
as many systems as are being created. 

The ASDs are only listed as possible because they all have q > 1 which, according to 
most prescriptions for mass transfer, will lead to common-envelope evolution and not steady 
transfer of material onto the WD (see Section 2). For these systems to become super-soft 
sources (and possibly SNIa) the common-envelope phase must be avoided, otherwise a DWD, 
or possibly a single giant, will result. The number of possible ASDs with q < 1 at the onset 
of mass transfer predicted by population synthesis is only ~ 0.4 per simulation which agrees 
with what we found. However, there is a fair amount of uncertainty in the value of g cr i t , 
the mass-ratio above which mass transfer from a giant proceeds on a dynamical timescale. 
Webbink (1988) presents an expression for g crit which differs by more than a factor of two 
compared to the model of Hurley, Tout & Pols (2002) for certain values of the giant envelope 
mass, and which gives q cr i t > 2 as the envelope mass becomes small. Furthermore, Yungelson 
& Livio (1998) have shown that the assumption of a strong optically thick wind from the 
accreting star can act to stabilize the mass transfer. If we allow all values of g cr it to lead to 
stable mass transfer (on a thermal timescale), and assume that all of the transferred material 
is accreted by the WD, then Table 3 shows that ~ 8 WDs per simulation will reach the Mch 
and explode. 
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5. 



Discussion 



The key result of this paper is that each of the possible progenitors of SNIa is pref- 
erentially manufactured in star clusters by orbital hardening and/or exchange interactions. 
Consequently, in the dense environment of a star cluster, it is binaries altered by, or created 
from, dynamical interactions that provide the dominant formation channel for short-period 
DWDs and possible ASDs. 



An overly simple prediction based on our results is that we might thus expect to see 
many or most Type la supernovae erupting in star clusters, often in their cores. There are 
two reasons to be cautious with this prediction. 

First, open star clusters disperse, usually on a 1 —6 Gyr timescale. Stars evaporate from 
clusters as they acquire energy from other cluster stars, as the Galactic tidal field incessantly 
tugs at cluster outliers, and through encounters with Giant Molecular Clouds. Second, we 
don't know what fraction of stars are created in star clusters, and what fraction of those 
escape during the earliest stages of cluster life. Observationally checking this prediction 
is difficult because only a dozen or so SNela have been identified closer than the Coma 
Cluster. Crowding will make at least some SNela appear "close" to star clusters, even if 
they are separated by hundreds of parsecs. 



At first glance, performing simulations of globular clusters should act to amplify the 
effects of dynamical interactions on the stellar populations within the cluster. These sim- 
ulations will operate at higher particle density and thus the rate of stellar encounters will 
increase. However, in some cases, this may have the effect of closing production channels. 
Consider the ASD systems which comprise a Roche-lobe filling giant star and a WD sec- 
ondary. For such a system to evolve to become a SNIa mass transfer onto the WD must 
be stable for a significant length of time. Observations of globular clusters (Guhathakurta 
et al. 1998, for example) show that bright giants are depleted in the core of the cluster 
indicating that they have been involved in collisions, a consequence of their relatively large 
cross-section. Such collisions could act to reduce the incidence of stable mass transfer from 
giants in a globular cluster, either by interrupting the mass transfer or preventing it from 
occurring at all. On the other hand, short-period DWD systems present a relatively small 
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cross-section for collision so the increased stellar density should lead to an increased num- 
ber of weak encounters. The orbital perturbations resulting from these will enhance the 
likelihood of DWD merger events. 

5.3. SNIa Birthrates 

The predicted Galactic birthrate of SNela from merging supra- Chandrasekhar DWDs 
is found by Hurley, Tout & Pols (2002) to be 2.6 x lO^yr -1 . This number is calculated 
from population synthesis of isolated binaries using the same input parameters to the binary 
evolution model as used in our NB0DY4 simulations. Alternatively, if the Gaussian period 
distribution for nearby solar-like stars found by Duquennoy & Mayor (1991) is used to 
determine the initial orbits of the binaries, this number changes by less than 4%. A sim- 
ilar rate is found by Tutukov & Yungelson (1994) using an independent binary evolution 
algorithm. The observed rate is 4 ± 1 x 10 _3 yr _1 (Cappellaro et al. 1997) so clearly a 
factor of 10 increase in the predicted rate would bring binary evolution models into conflict 
with observations. However, Hurley, Tout & Pols (2002) found that if they draw the binary 
component masses independently from a single star IMF, rather than using a uniform mass- 
ratio distribution, the predicted rate drops by at least a factor of 10. So the enhancement 
of SNIa candidates found in our open cluster simulations would allow agreement for the 
more general mass-ratio distribution when trying to match the results of binary population 
synthesis studies to observations. 

There are additional parameters intrinsic to the binary evolution algorithm which can 
also affect the production rate of DWDs. The common-envelope efficiency parameter is a 
prime example. Reducing q:ce from 3.0 to 1.0 makes it harder to form short-period binaries 
via the common-envelope channel and indeed, Hurley, Tout & Pols (2002) found that this 
also decreased the predicted birthrate of SNela from merging supra- Chandrasekhar DWDs 
by a factor of 10. Because the critical mass-ratio for dynamical mass transfer affects the 
frequency of common-envelope events, uncertainty in its value (see previous section) will 
similarly affect any predicted number. Trying to constrain the input parameters to a model 
of binary evolution by matching the results of binary population synthesis to observations 
of particular stellar populations is a risky business. The number of uncertain parameters is 
large and it is entirely possible that an error in one value will mask the error in another. 
But the key result of this paper remains: hardening of binaries occurs only in clusters, and 
this preferentially creates SNIa candidates. 
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5.4. Binary Fraction 

Clearly the number of DWD binaries produced in the simulations would increase if we 
had used a higher binary fraction. The fraction of 10% is a rather conservative choice and 
we note that in some open clusters it has been found that binaries may be as, or even more, 
populous than single stars (Fan et al. 1996, for example). We cannot stress enough however 
(this goes for the discussion in the previous paragraph as well) that it is the number of 
close DWDs produced in the simulations relative to the number produced in the field that 
is of primary interest in this work. Binaries do present a greater cross-section for dynamical 
interaction than do single stars (Heggie, Hut & McMillan 1996) so a higher binary fraction 
does have the capacity to affect the simulation results. This will be addressed in future work. 

5.5. Eccentricity 

The treatment of gravitational radiation used in the binary evolution algorithm contains 
a strong dependence on orbital eccentricity, i.e. the removal of angular momentum from the 
short-period system is accelerated if the orbit is eccentric. Binaries emerging from a common- 
envelope phase are assumed to be circular so the only way for a short-period DWD orbit 
to be eccentric is through dynamical interactions subsequent to DWD formation. However, 
inducing an eccentricity in to an already circular orbit, especially in the case of small orbital 
separation, is difficult (Heggie & Rasio 1996). In none of the four simulated SNIa candidates 
that required perturbations to the orbit of a DWD was a detectable eccentricity induced as 
part of the hardening process. 

On a related matter, Hurley et al. (2001) have shown that the shape of the eccentricity 
distribution used for the primordial binary population can affect the number of blue stragglers 
produced from that population, for example. This is less important in the case of DWD 
binaries as the systems of interest will simply form from a slightly different set of primordial 
binaries. The same goes for variations in the way that tidal friction is modelled (see Hurley, 
Tout & Pols 2002 for a detailed discussion). 

5.6. DWDs and CVs in Clusters 

Another intriguing suggestion from this work is that DWDs of all periods should be 
rather commonplace in open clusters, and likely segregated towards the cluster centres. A 
significant fraction of these DWDs must have been involved in exchange interactions. We 
encourage observers to survey for such objects. In the case of CVs we do not find a similar 
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enhancement in open clusters: typically 1 — 2 CVs were formed in each simulation. We 
expect that this is mainly due to the required presence of a low-mass MS star, less massive 
than its WD companion, in any binary that will evolve to a long-lived, and thus stable, 
CV state. This means that the progenitors of CVs will on average be less massive than 
the progenitors of DWDs, either on the ZAMS or after the WD (or WDs) have formed, so 
that they are less likely to reside in the core of the cluster and be involved in dynamical 
interactions. Furthermore, in a 3-body exchange interaction, it is normally the least massive 
star that is ejected so CVs are unlikely to form in this manner. We do not rule out the 
possiblity that CV production will be enhanced in simulations of globular clusters. We also 
note that if merging supra-Chandrasekhar DWDs lead to the formation of neutron stars via 
an AIC then they are still extremely interesting objects in terms of the neutron star retention 
problem in star clusters (Pfahl, Rappaport & Podsiadlowski 2002). 



6. Conclusions 

We find a greatly increased rate of production of "loaded guns" - supra-Chandrasekhar 
double-white-dwarfs with inspiral ages shorter than a Hubble time - in star clusters relative 
to the field. Orbital hardening and exchange interactions are the responsible mechanisms 
for the enhanced rates. Neither of these processes operates in the field, and neither has been 
included in previous population studies of SNIa progenitors. 

The production rate of possible accreting single degenerates is not significantly enhanced 
relative to the field. However, a major fraction of these systems are formed in exchange 
interactions, which means that modelling of the stellar dynamics is destroying production 
channels as well as creating them. Whether or not accreting single degenerates will evolve to 
become super-soft sources, and ultimately type la supernovae, depends critically on avoiding 
the phase of common-envelope evolution that occurs if mass transfer proceeds on a dynamical 
timescale. 

The results are based on studies of open cluster size iV-body simulations, but we expect 
the effect to be even stronger in the case of globular clusters. The more frequent encounters 
in globulars will harden DWDs much more rapidly than in open clusters, and hence increase 
the predicted rate of creation of SNIa progenitors. An increased rate of collisions involving 
giant stars will deplete the numbers of accreting single degenerate binaries. 

We are extremely grateful to Jun Makino and the University of Tokyo for the loan of 
the GRAPE-6 board. The generous support of Doug Ellis and the Cordelia Corporation has 
enabled AMNH to purchase new GRAPE-6 boards and we are most grateful for that. MS 



-14- 



thanks Ken Nomoto for very helpful discussions and suggestions. We thank Fred Rasio and 
Vicky Kalogera for organizing the 2001 Aspen Center for Physics meeting on Star Clusters 
where some of the ideas for this paper were clarified. We also thank the referee for comments 
that greatly improved certain aspects of this manuscript. 



- 15 - 



REFERENCES 

Aarseth, S., Henon, M., & Wielen, R. 1974, A&A, 37, 183 
Aarseth, S. J. 1999, PASP, 111, 1333 

Cappellaro, E., Turatto, M., Tsvetkov, D.Yu., Bartunov, O.S., Pollas, C, Evans, R., & 
Hamuy, M. 1997, A&A, 322, 431 

Chen, K, & Leonard, P.J.T. 1993, ApJ, 411, L75 

Chernoff, D.F., & Weinberg, M.D. 1990, ApJ, 351, 121 

Duquennoy, A., & Mayor, M. 1991, A&A, 248, 485 

Eggleton, P.P. 2002, Evolutionary Processes in Binary and Multiple Stars, (Cambridge: 
Cambridge University Press), in preparation 

Eggleton, P.P., Fitchett, M., & Tout C.A. 1989, ApJ, 347, 998 

Fan, X., et al. 1996, AJ, 112, 628 

Giersz, M., & Heggie, D.C. 1997, MNRAS, 286, 709 

Guhathakurta, P., Webster, Z.T., Yanny, B., Schneider, D.P, & Bahcall, J.N. 1998, AJ, 116, 
1757 

Heggie, D.C. 1975, MNRAS, 173, 729 

Heggie, D.C, Hut, P., & McMillan, S.L.W. 1996, ApJ, 467, 359 

Heggie, D.C, & Rasio, F.A. 1996, MNRAS, 282, 1064 

Hurley, J.R., Pols, O.R., & Tout, C.A. 2000, MNRAS, 315, 543 

Hurley, J. R., Tout, C. A., Aarseth, S. J., & Pols, O.R. 2001, MNRAS, 323, 630 

Hurley, J.R., Tout, C.A., & Pols, O.R. 2002, MNRAS, 329, 897 

Hurley, J.R., & Shara, M.M. 2002, ApJ, in press 

Iben, I.- Jr., & Livio, M. 1993, PASP, 105, 1373 

Kraft, R. P. 1983, Highlights in Astronomy, 6, 129 

Kroupa, P., Tout, C. A., & Gilmore, G. 1991, MNRAS, 251, 293 



-16- 



Kroupa, P., Tout, C. A., & Gilmore, G. 1993, MNRAS, 262, 545 

Kurucz R.L. 1992, in IAU Symposium 149, The Stellar Populations of Galaxies, ed. B. 
Barbuy, & A. Renzini (Dordrecht: Kluwer), 225 

Lada, E. A., Strom, K. M., & Myers, P. C. 1993, in Protostars and Planets III, ed. E. Levy, 
& J. Lunine (University of Arizona), 245 

Leibundgut, B. 2001, ARA&A, 39, 67 

Makarov, V.V., & Fabricius, C. 2001, A&A, 368, 866 

Makino, J. 2001, in ASP Conf. Ser. XX, Stellar Collisions, Mergers and their Consequences, 
ed. M. M. Shara (San Francisco: ASP), in press 

Mardling, R.A., & Aarseth, S.J. 2001, MNRAS, 321, 398 

Perlmutter, S., et al. 1999, ApJ, 517, 565 

Pfahl, E., Rappaport, S., & Podsiadlowski, Ph. 2002, ApJ, submitted 

Phillips, M.M. 1993, ApJ, 413L, 105 

Portegies Zwart, S.F., & Verbunt F. 1996, A&A, 309, 179 

Riess, A.G. et al. 1998, AJ, 116, 1009 

Riess, A.G. 2000, PASP, 112, 1284 

Saffer, R.A., Livio, M., & Yungelson, L.R. 1998, ApJ, 502, 394 
Saio, H., & Nomoto, K. 1998, ApJ, 500, 388 

Tout, C.A., Aarseth, S.J., Pols, O.R., & Eggleton P.P. 1997, MNRAS, 291, 732 

Tout, C.A., Regos, E., Wickramasinghe, D., Hurley, J.R., & Pols, O.R. 2001, in ASP Conf. 
Ser. 229, Evolution of Binary and Multiple Star Systems: A Meeting in Celebration 
of Peter Eggleton's 60th Birthday, ed. Ph. Podsiadlowski, S. Rappaport, A. R. King, 
F. D'Antona, & L. Burder, (San Francisco: ASP), 275 

Tutukov, A., & Yungelson, L. 1994, MNRAS, 268, 871 

Tutukov, A., & Yungelson, L. 1994, MNRAS, 280, 1035 

Webbink R.F. 1998, in IAU Colloquim 103, The Symbiotic Phenomenon, ed. J. Mikolajewska, 
M. Friedjung, S.J. Kenyon, & R. Viotto (Dordrecht: Kluwer), 311 



-17- 

Yungelson, L., & Livio, M. 1998, ApJ, 497, 168 
Yungelson, L., & Livio, M. 2000, ApJ, 528, 108 

Yungelson, L., Livio, M., Truran, J.W., Tutukov, A., & Federova, A.V. 1996, ApJ, 466, 890 



This preprint was prepared with the AAS IATgX macros v5.0. 



-18- 



Table 1. Double- WD systems that are Type la candidates. To qualify the system must 
have a combined mass in excess of the Chandrasekhar mass, 1.44M , and a gravitational 
radiation merger timescale less than the age of the Universe, ~ 1.2 x 10 10 yr. The time at 
formation for the double- WD (DWD) system, Myr units, is given in Column I. The types 

of the WDs are listed in Column 2. Three types of WD are distinguished: helium 
composition (He), carbon-oxygen (CO), and oxygen-neon (ONe). The individual masses of 
the two WDs are given in Columns 3 and 4 respectively, and the combined mass is given in 
Column 5. All masses are in solar units. The period of the binary is given in Column 6 in 

units of days. Column 7 gives an estimate of the time it will take the DWD system to 
merge, in yrs, owing to angular momentum loss from gravitational radiation. This estimate 
comes from integrating eq. (48) of Hurley, Tout & Pols (2002). Column 8 summarizes the 
history of each system using the following code: primordial binary (PRIM); exchange 

interaction (EXCH); perturbation before (PB), or after (PA), Double- WD formed; 
subsequent escape from cluster (ESC); subsequent disruption in dynamical encounter 
(DISR). Note that perturbations to the orbit are only recorded if they lead to a change in 

the evolution path of the binary. 
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Table 2. Double- WD systems that are not Type la candidates. Formed either from a 
primordial binary (PRIM) or via an exchange (EXCH) interaction. See Table 1 for an 

explanation of what each column entails. 



Tform Types M x M 2 M h Period T grav Legend 
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Table 2 — Continued 
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Table 2 — Continued 
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Table 3. Potential super-soft sources. Time at onset of mass transfer is given in 
Column 1, in Myr units, and the stellar types of the stars at this time are listed in 
Column 2: Roche-lobe filling star is either on the Hertzsprung gap (HG) or the first giant 
branch (GB). The mass of the sub-giant or giant, the mass of the WD, and the mass-ratio 
(q = Mi/M 2 )., are given in Columns 3, 4, and 5, respectively. Column 6 shows the resultant 
mass of the WD if it were to accept all of the mass available in the donor stars envelope at 
the onset of mass-transfer. All masses are in solar units. The period of the binary is given 
in Column 7 in units of days. Column 8 summarizes the history of each system using the 
following code: primordial binary (PRIM); exchange interaction (EXCH); perturbation to 

orbit before Roche- lobe overflow (PB); steady mass-transfer (SMT); common-envelope 
evolution (CE). The outcome of CE is the creation of a DWD system, where the giant has 
become a He WD, but in cases where the giant core was non-degenerate a naked helium 

star (HeMS) is the resulting companion. 
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Fig. 1. — The evolution of particular binary systems plotted as the logarithm of the absolute 
value of the binding energy, E\>, where E\> is in units of Mq/AU. The top panel shows 
the ONe-CO DWD formed at 334 Myr in the first Z = 0.004 simulation. This system is a 
possible SNIa candidate (fourth entry in Table 1) and is a primordial binary that had its 
orbit perturbed after the DWD formed. The system ceases to exist when, at 630 Myr, the 
two WDs merge to form a supra- Chandrasekhar object. The lower panel shows the ONe-CO 
DWD formed at 1 192 Myr in the second Z = 0.02 simulation. This system is also a possible 
SNIa candidate (last entry in Table 1) and formed in an exchange interaction at 620 Myr. 
The evolution of the primordial binary involved in the exchange is shown as a dashed line. 
The evolution after 3 000 Myr, in which the two WDs continue to spiral together owing to 
gravitational radiation, is omitted for clarity. 
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Fig. 2. — Histogram of double- WD 
shown. All systems listed in Tables 1 



masses. The Chandrasekhar mass of 1.44M is also 
and 2 are included. 
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Fig. 3. — Histogram of double- WD periods at the time of formation. All systems listed in 
Tables 1 and 2 are included. 
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Fig. 4. — Colour- magnitude diagrams for the Z = 0.004 simulations at 1.0, 2.0, 3.0 and 
4.0 Gyr. Stars from both simulations are plotted but the numbers of single stars and bina- 
ries given in each figure are averaged. Main-sequence stars (dots), blue stragglers (stars), 
sub-giants, giants and naked helium stars (open circles) and white dwarfs (dots) are distin- 
guished. Binary stars are denoted by overlapping symbols appropriate to the stellar type 
of the components, with main-sequence binary components depicted with filled circles and 
white dwarf binary components as diamonds. Type la candidates are circled. Bolometric 
corrections computed by Kurucz (1992) from synthetic stellar spectra are used to convert 
theoretical stellar quantities to observed colours. These corrections are strictly not valid for 
WDs and extremely cool giants. 
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Fig. 5. — Same as Figure 4 but for the Z = 0.02 simulations. Note that many of the single 
main-sequence stars are overlayed by binary stars (in Figure 4 as well). 
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Fig. 6. — Population gradients in the XY-plane for single stars, binaries and double- WDs, 
averaged over the four simulations presented in this paper. The tidal radius of the cluster is 
shown as a vertical dashed line. Note that stars are not actually removed from a simulation 
until they are at a distance greater than two tidal radii from the cluster centre. The half-mass 
radius, r h , the tidal radius, and the numbers of each sub-population, are averaged values per 
simulation. 



